Geog533 Lab 10

Question 1

Use R package UScensus2010county to complete the following tasks: (20 pt.)

if(!require(UScensus2010)) {install.packages("UScensus2010")}

Loading required package: UScensus2010
Loading required package: maptools
Checking rgeos availability: TRUE
Loading required package: foreign


Package UScensus2010: US Census 2010 Suite of R Packages
Version 0.11 created on 2011-11-18.

Zack Almquist, University of California-Irvine
ne

For citation information, type citation("UScensus2010").
Type help(package=UScensus2010) to get started.

library(UScensus2010)
if(!require(UScensus2010county)) {install.county("osx")}

Loading required package: UScensus2010county

UScensus2010county: US Census 2010 County Level Shapefiles and Additional
Demographic Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package="UScensus2010county") to get started.

For citation information, type citation("UScensus2010county").

library(UScensus2010county)
library(tmap)

package ‘tmap’ was built under R version 3.4.3

Question 1(a)

Plot a map of New York counties using the plot function.

data(new_york.county10)
shpfile <- new_york.county10
df1 <- shpfile@data
df1

plot(shpfile)

Question 1(b)

Plot a map of New York counties using the qtm function.

install.packages("tmap")

Error in install.packages : Updating loaded packages

library(tmap)
qtm(shpfile)

Question 1(c)

How many counties in New York State?

df

function (x, df1, df2, ncp, log = FALSE) 
{
    if (missing(ncp)) 
        .Call(C_df, x, df1, df2, log)
    else .Call(C_dnf, x, df1, df2, ncp, log)
}
<bytecode: 0x1074b8040>
<environment: namespace:stats>

print("There are 62 counties in New York State")

[1] "There are 62 counties in New York State"

Question 1(d)

What’s the 3-digit fips code of Broome County?

print("36007")

[1] "36007"

Question 1(e)

Compute descriptive statistics of the population column (P0010001), including total, minimum, maximum, mean, median, and skewness.

library(moments)
nypop <- df$P0010001

Error in df$P0010001 : object of type 'closure' is not subsettable

Question 1(f)

Create a histogram and a boxplot of the population.

rr a <- hist(nypop)

Error in hist.default(nypop) : 'x' must be numeric

Question 2

Use R package UScensus2010tract to complete the following tasks: (20 pt.)

rr if(!require(UScensus2010tract)) install.tract()

Loading required package: UScensus2010tract
there is no package called ‘UScensus2010tract’trying URL 'http://lakshmi.calit2.uci.edu/census2000/R/src/contrib/UScensus2010tract_1.00.tar.gz'
Content type 'application/x-gzip' length 190473496 bytes (181.6 MB)
==================================================
downloaded 181.6 MB

* installing *source* package ‘UScensus2010tract’ ...
** R
** data
** inst
** preparing package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded

UScensus2010tract: US Census 2010 Tract Level Shapefiles and Additional Demographic
Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package=\UScensus2010tract\) to get started.

For citation information, type citation(\UScensus2010tract\).
* DONE (UScensus2010tract)

The downloaded source packages are in
    ‘/private/var/folders/l3/pkmg9p591775dgyt86k3s9ym0000gn/T/RtmpIedmm3/downloaded_packages’
NULL

rr library(UScensus2010tract)

UScensus2010tract: US Census 2010 Tract Level Shapefiles and Additional Demographic
Data 
Version 1.00 created on 2011-11-06 
copyright (c) 2011, Zack W. Almquist, University of California-Irvine
Type help(package=\UScensus2010tract\) to get started.

For citation information, type citation(\UScensus2010tract\).

Question 2(a)

Plot a map of New York census tracts using the plot function.

rr data(new_york.tract10) shpfile <- new_york.tract10 df <- shpfile@data df r plot(shpfile)

Question 2(b)

Compute the total population based on census tracts.

df

function (x, df1, df2, ncp, log = FALSE) 
{
    if (missing(ncp)) 
        .Call(C_df, x, df1, df2, log)
    else .Call(C_dnf, x, df1, df2, ncp, log)
}
<bytecode: 0x10d7b5578>
<environment: namespace:stats>

newpop <- df$P0010001

Error in df$P0010001 : object of type 'closure' is not subsettable

Question 2(c)

Select all census tracts in Broome County and plot the map.

broome <- county(fips="007",name="Broome",state="NY",level="tract")

Error in county(fips = "007", name = "Broome", state = "NY", level = "tract") : 
  could not find function "county"

Question 2(d)

What’s the total population of Broome County?

rr bpop <- sum(broome$P0010001) bpop

[1] 200600

Question 2(e)

Create a histogram and a boxplot of population based on census tracts of Broome County.

rr hist(broome$P0010001)

rr boxplot(broome$P0010001)

Question 2(f)

Select the first five columns of the shapefile of Broome County census tract; add a new population ratio column (= census tract population / county population); save the new shapefile to the hard drive.

popra <- data.frame(broome[,1:5])

Error in data.frame(broome[, 1:5]) : object 'broome' not found

Question 3

Use R packages raster and ncdf4 to complete the following tasks: (20 pt.)

library(raster)
library(ncdf4)

Question 3(a)

Load the multi-band raster dataset “NDVI.nc” into R.

ndvi <- raster("NDVI.nc")

Question 3(b)

Get the basic information about the dataset, including the number of rows, columns, cells, and bands; spatial resolution, extent, bounding box, and projection.

nrow(ndvi)

[1] 360

ncol(ndvi)

[1] 720

ncell(ndvi)

[1] 259200

nbands(ndvi)

[1] 36

res(ndvi)

[1] 0.5 0.5

extent(ndvi)

class       : Extent 
xmin        : -180 
xmax        : 180 
ymin        : -90 
ymax        : 90 

bbox(ndvi)

    min max
s1 -180 180
s2  -90  90

projection(ndvi)

[1] "+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0"

Question 3(c)

Aggregate all bands to generate a mean NDVI raster and a maximum NDVI raster; save the two new raster datasets to the hard drive.

ndvi.rb <- brick("ndvi.nc")
meanndvi <- mean(ndvi.rb)
plot(meanndvi)

writeRaster(meanndvi,filename = "meanndvi.tif",overwrite=TRUE)
maxndvi <- max(ndvi.rb)
plot(maxndvi)

writeRaster(maxndvi,filename = "maxndvi.tif",overwrite=TRUE)

Question 3(d)

Plot the maps of monthly NDVI of the year 2001

ndvi2k <- ndvi.rb[[13:24]]
plot(ndvi2k)

Question 3(e)

Create histograms of monthly NDVI of the year 2001.

hist(ndvi2k)

size changed to n because it cannot be larger than n when replace is FALSE39% of the raster cells were used (of which 42% were NA). 58093 values used.

Question 3(f)

Plot the NDVI map of July 2000; click any location with data on the map and retrieve the NDVI time series for all years; plot the NDVI time series of the selected location.

plot(ndvi.rb,7) 
location <- ndvi.rb[150, 50]
location <- click(ndvi.rb, n=1, xy=TRUE)

location <- click(ndvi.rb, n=1, xy=FALSE)